
import matplotlib.pyplot as plt
import cobra
from cobra.io import validate_sbml_model
import importlib
from xml.etree import ElementTree
import utils.Model_correction as mc
import sys
import utils.model_maj as mj
import utils.viz_utils as vu
#import cplex
import cbmpy
import pandas as pd
import plotly.express as px
from plotly.subplots import make_subplots
import plotly.graph_objects as go
from tqdm import tqdm
from itertools import cycle
import seaborn as sns

pyparsing import

INFO: No xlwt module available, Excel spreadsheet creation disabled
CBGLPK based on swiglpk: not all methods implimented yet! 5.0

*****
Using CPLEX
*****

doFBAMinSum not available with GLPK

INFO: No xlrd module available, Excel spreadsheet reading disabled


***********************************************************************
* Welcome to CBMPy (0.8.4) - PySCeS Constraint Based Modelling        *
*                http://cbmpy.sourceforge.net                         *
* Copyright(C) Brett G. Olivier 2014 - 2020                           *
* Systems Biology Lab, Vrije Universiteit Amsterdam                   *
* Amsterdam, The Netherlands                                          *
* CBMPy is developed as part of the BeBasic MetaToolKit Project       *
* Distributed under the GNU GPL v 3.0 licence, see                    *
* LICENCE (supplied with this release) for details                    *
***********************************************************************


HepG2, errors = validate_sbml_model("../models_storage/Hep_G2_v7.xml")


iHep, errors = validate_sbml_model("../models_storage/iHep_v3.xml")


HepG2_C = HepG2.copy()


HepG2.reactions.EX_m01965x.bounds = (-565.0,-565.0) #Glucose exchange, set to the maximal input value observed with FVA
#iHep.reactions.EX_m01965x.bounds = (-1000.0, -1000.0)
iHep.reactions.HMR_4281.bounds = (0.0,0.0)


HepG2.reactions.HMR_3883.bounds = (0.0,0.0) #Serine --> Pyruvate exchange
HepG2.reactions.EX_m02630x.bounds = (0.0,0.0) #O2 exchange
HepG2.reactions.EX_m02819x.bounds = (0.0,1000.0) #Pyruvate exchange
HepG2.reactions.EX_m01910x.bounds = (0.0,0.0) #Galactose exchange
HepG2.reactions.EX_m01743x.bounds = (0.0, 1000.0) #D-Ribulose exchange
HepG2.reactions.EX_m01962x.bounds = (0.0,1000.0) #glucosamine exchange
HepG2.reactions.HMR_4316.bounds = (-1000.0,1000.0) # Glucose --> D-Glucitol
HepG2.reactions.EX_m02896x.bounds = (0.0,1000.0) # Serine intake
HepG2.reactions.EX_m01682x.bounds = (-1000.0,0.0) # Glucitol secretion blocked. Kept the intake just in case.
HepG2.reactions.EX_m01840x.bounds = (-1000.0,0.0) #Fructose exchange
#HepG2.reactions.HMR_4930.bounds = (-1000.0,1000.0) # Pyruvate transfer from cytoplasm to peroxysome
#HepG2.reactions.HMR_4281.bounds = (0.0,0.0) # Fermentation in peroxysome


iHep.objective = "biomass_components"
sol_iHep = iHep.optimize()
HepG2.objective = "biomass_components"
sol_G2 = HepG2.optimize()
HepG2_C.objective = "biomass_components"
sol_G2_C= HepG2_C.optimize()


df_Ci, df_Ri, df_Si, df_Mi, df_Pi, df_Xi, df_Li, df_Gi, df_Ni = vu.build_reaction_df(iHep)
df_C, df_R, df_S, df_M, df_P, df_X, df_L, df_G, df_N = vu.build_reaction_df(HepG2)
subS_HepG2 = vu.get_subsystem_fluxes(vu.build_reaction_df(HepG2))
subS_iHep = vu.get_subsystem_fluxes(vu.build_reaction_df(iHep))
subS_G2C = vu.get_subsystem_fluxes(vu.build_reaction_df(HepG2_C))
df_both = pd.concat([subS_HepG2, subS_G2C, subS_iHep], axis = 1)
df_both.columns = ["HepG2","HepG2_copy", "iHep"]
df_both = df_both.loc[(df_both["iHep"] != 0.0)& (df_both["HepG2"] != 0.0)]

clustermap = sns.clustermap(df_both.fillna(0), figsize =(10,20))


barplots_dict = vu.compartment_fluxes_barplots(iHep, HepG2)


title = "Cytoplasm reactions -- <b>CANCER MODEL<br />"
vu.plot_treemap(df_C,HepG2,title).show(renderer = 'notebook')


title = "Cytoplasm reactions -- <b>HEALTHY MODEL <br />"
vu.plot_treemap(df_Ci, iHep, title).show(renderer = "notebook")


barplots_dict["C_c"][0]


title="Mitochondrial reactions -- <b>HEALTHY MODEL<br />"
vu.plot_treemap(df_Mi, iHep, title).show(renderer = "notebook")


title="Mitochondrial reactions -- <b>CANCER MODEL<br />"
vu.plot_treemap(df_M, HepG2, title).show(renderer = "notebook")


barplots_dict["C_m"][0]


title="Peroxysomal reactions -- <b>HEALTHY MODEL<br />"
vu.plot_treemap(df_Pi, iHep, title).show(renderer = "notebook")


title="Peroxysomal reactions -- <b>CANCER MODEL<br />"
vu.plot_treemap(df_P, HepG2, title).show(renderer = "notebook")


barplots_dict["C_p"][0]

Models loading and bounds modification¶

Optimisation¶

Heatmap¶

Comparaison, par compartiment, des subsystems les plus actifs pour les deux modèles :¶

Treemaps + barplots¶

Cytoplasme¶

Mitochondrie¶

Peroxysome¶